home *** CD-ROM | disk | FTP | other *** search
/ Language/OS - Multiplatform Resource Library / LANGUAGE OS.iso / dsp / fft / fft_eyal.lha / fft_eyal / readme < prev   
Encoding:
Text File  |  1991-09-01  |  6.3 KB  |  141 lines
  1. This is the second release of these programs. The first release (packaged
  2. by Udi Finkelstein) included only the m68k programs.
  3.  
  4. This is my implementation of integer real fft computation. It was
  5. designed for high speed analysis of audio input in real time on a
  6. PC/AT. It has support for generic C and optimized assembly for:
  7. - intel 80x86
  8. - motorola 68k
  9. - natsemi 32k
  10. The assembler is usualy 4-8 times faster than the generic C so it is worth
  11. using.
  12.  
  13. The data used was 8bit but the routine does 16bit signed arithmetic.
  14. The algorithm of choice was one which minimizes mults while keeping the
  15. number of needed registers to a minimum (the Intel chip does not have
  16. enough). Data if gradualy scaled in flight to avoid overflow.
  17.  
  18. A standard (2,4) split-radix was used as described in [Sorensen et al.,
  19. IEEE ASSP-35 no 6., June 87, pp. 849-863]. Other algoriths (WFT etc.)
  20. where found to use too complex arithmetic for efficient implementation
  21. on a machine with a decent mult instruction because register spills slow
  22. it more than the reduce mults speed it.
  23.  
  24. WHAT IT DOES NOT DO:
  25.  
  26. - NO inverse fft.
  27. - NO complex input.
  28. - NO windowing.
  29. - output is power spectrum WITHOUT the square root taken.
  30. - the sample data array (x[]) is not descramled so you cannot use it directly.
  31.   see how the power spectrum is derived if you want to generate phase or
  32.   whatever. fft7() and fft8() function do this post-processing.
  33.  
  34. The method of implementaion was to have a program (fftg.c) that
  35. generates the fft routine which is then compiled/assembled and used.
  36. The generated program is assembly, but for portability this release
  37. generates a C program too. I have routines for Intel assembly and ns32k
  38. assembly, other machines can have the basic routines re-written without
  39. much effort (but for best performance one would optimize these routines
  40. extensively).
  41.  
  42. The archive contains a sample program (fft2.c) that uses the fft
  43. routine. It draws some basic picture on a screen using move() and
  44. draw() function which you will have to supply (it uses the microsoft c
  45. graphics library in hercules mode using msherc.com). To run the test
  46. program: 'fft2a' or 'fft2c' or 'fft2f'. you can try 'fft2a x' which
  47. will call fft 10000 times and display the time, then do 'fft2a x x' to
  48. measure the overhead of the loop, then a result of (e.g.) 29 from the
  49. first and 3 from the second means an average time of 2.6ms/fft. When
  50. the program asks for a data file respond 'd1' of 'd2' etc. or 'end' to
  51. stop.
  52.  
  53. Originaly the program was a BASIC program from some British mag a few
  54. years ago which started my interest; it took a few minutes to do one
  55. fft. The current version does a 256 point fft on a 25MHz ns32532 in
  56. 2.6ms(asm)/4.2ms(c)/5.8ms(funcs), on a 10MHz 80286 it takes 7ms(asm, it
  57. has a rather fast int mult)/38.8ms(c)/44.4ms(funcs). As you see, gcc on
  58. the ns32k compares better to hadcrafted assembly than msc5.1.
  59.  
  60. I now use a 20MHz 80286. If you sample at 44Khz, 256 points take 5.8ms
  61. to sample but only 2.9ms to do the fft (assuming the sampling is done
  62. via DMA, inerrupts would eat too much cpu time), so you have 2.9ms for
  63. other stuff (display? this is slow unless you really try. I wrote my
  64. low level display routines or the fft gain would be lost). On a 386/486
  65. you are laughing.
  66.  
  67. The package contains one ready fft() in the file fft8c.f; it does a 256
  68. point fft in the slowest way but it should compile readily.
  69.  
  70. A program to measure the time is also supplied as fft1.c. You will need to
  71. link the following:
  72.     fft1.o fft8f.o fftsubs.o
  73. then run the result to see how fast it goes. Be patient, the test takes about
  74. one minute. On some machine (where 'mult' has variable time) this test will
  75. give only the typical speed.
  76.  
  77. To get the generator programs 'make progs'. A generator program is built
  78. from two modules: fft.g (the core) and fftout?.c (the output driver). You
  79. should choose one of the drivers or write your own if there is no support
  80. for your machine (in which case I would like to hear from you). The driver
  81. fftoutf.c generates a C program that calls functions in fftsubs.c; this
  82. should work on every system.
  83.  
  84. Once you have the generator, you need to run it and generate the fft program:
  85.     fftgf fft8f 8 fft
  86. this will build a 256 (2^8) points fft() function into fft8f.c. Now compile
  87. and then link into your application.
  88.  
  89. All of the above is done in the makefile, look at it.
  90.  
  91. To get the fft programs 'make fft'. The 'c' version uses fftsubs.h. The
  92. 'f' version uses C function calls which is a bit slower but uses far
  93. less memory; it uses fftsubs.c. The '86' version uses fftsubs.mac, Note
  94. that the makefile generates a 256 point routine with entry point fft(),
  95. you may want to change it for your needs or just edit the output.
  96.  
  97. fftouts.c generates ns32k assembly, fftouta.c generates Intel assembly.
  98. These should give a good idea of how to port the assembly level
  99. routines to another machine.
  100.  
  101. To get the test programs 'make tests' (but first fix fft2.c for the
  102. graphics if needed).
  103.  
  104. MSDOS NOTE: The c output has to be split into smaller files because a)
  105. the compiler wont optimize a large prog b) the compiler dies on any
  106. larger progs, and c) the code segment would be too large. The output
  107. routines in fftoutc.c will break the output every 1000 statements. This
  108. is fine for msc5.1, other compilers may cope with larger programs, so
  109. just up the limit in break_fft().
  110.  
  111. HOW TO MODIFY FFT:
  112.  
  113. The program to modify is fftg.c and the related fftout?.c that you use.
  114. You would probably modify the latter.
  115.  
  116. In general, the calls to fft1...fft5 handle the fft proper while
  117. fft7+fft8 do the final power spectrum and reordering. The last two is
  118. where you extract the results you wanted from x[] into qf[]. x[] holds
  119. the real and complex parts in the order:
  120.     xx[0-255]= R[0],...,R[128],I[127],...,I[1]
  121. but in this program xx[i] is in position x[mm[i]].
  122.  
  123.  
  124. HOW TO CALL FFT:
  125.  
  126. In your program header include:
  127.     short near x[256] = {0};    /* some compilers need to initialize */
  128.     short near qf[128+1] = {0};
  129.     extern far fft ();    /* or whatever you called it */
  130. Then in the program:
  131.     /* read the data points into x[0-255]*/
  132.     fft ();
  133.     /* the power spectrum is in qf[0-128], qf[0] is the DC componnent */
  134. Thats all.
  135.  
  136. Please let me know of problems/difficulties you encounter with these
  137. programs as well as any suggestions/improvements.
  138.  
  139. Regards
  140.     Eyal Lebedinsky        eyal@ise.canberra.edu.au
  141.